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Giant exoplanets orbiting close to their host stars are unlikely to have formed in 
their present configurations’. These ‘hot Jupiter’ planets are instead thought to have 
migrated inward from beyond the ice line and several viable migration channels 
have been proposed, including eccentricity excitation through angular-momentum 
exchange with a third body followed by tidally driven orbital circularization”*. The 
discovery of the extremely eccentric (e = 0.93) giant exoplanet HD 80606 b (ref. 4) 
provided observational evidence that hot Jupiters may have formed through 

this high-eccentricity tidal-migration pathway’. However, no similar hot-Jupiter 
progenitors have been found and simulations predict that one factor affecting the 
efficacy of this mechanism is exoplanet mass, as low-mass planets are more likely to 
be tidally disrupted during periastron passage? ®. Here we present spectroscopic and 
photometric observations of TIC 241249530 b, a high-mass, transiting warm Jupiter 
with an extreme orbital eccentricity of e = 0.94. The orbit of TIC 241249530 bis 
consistent with a history of eccentricity oscillations and a future tidal circularization 
trajectory. Our analysis of the mass and eccentricity distributions of the transiting- 
warm-Jupiter population further reveals a correlation between high mass and high 


eccentricity. 


The Transiting Exoplanet Survey Satellite (TESS)? monitored the 
apparent brightness of the star TIC 241249530 for 28 days during the 
second year of its primary mission. These data reveal a transit-like 
approximately 0.8% dip in brightness on 12 January 2020, the shape 
and depth of which were consistent with a Jupiter-sized planet pass- 
ing in front of the star (Fig. 1a). To find out the nature and origin of 
this signal, we conducted a series of ground-based observations 
of TIC 241249530. We first used high-spatial-resolution speckle 
imaging data from NESSI” to rule out the presence of contaminat- 
ing sources and confirm that the signal was not associated with a 
background eclipsing-binary in the TESS aperture. We then began 
radial velocity (RV) observations with the NEID spectrograph”, which 
revealed that the TESS transit was probably induced by a giant exo- 
planet on a highly eccentric (e = 0.94), long-period (P = 167 days) 
orbit. These measurements were consistent with the absence of a 


transit detection when TESS re-observed this star for 27 days from 
December 2022 to January 2023. Further NEID measurements, sup- 
plemented by observations with the HPF” and HARPS-N® spectro- 
graphs, were strategically scheduled to be taken when the planet 
was predicted to be approaching periastron and thus inducing large 
stellar RV variations. We attempted to detect a second transit using 
the global Unistellar telescope network” in March 2023, but these 
efforts were unsuccessful as the ephemeris was not yet well con- 
strained. However, RV data collected during the periastron window 
enabled us to more precisely predict the subsequent transit win- 
dow. We captured the first half of this transit using the engineered 
diffuser” on the ARCTIC imager” on 30 August 2023 (Fig. 1b). We 
refined the orbit using the ARCTIC data together with further NEID 
observations, including several concurrent with this transit and the 
subsequent one on 12 February 2024. Our ensemble of photometric 
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Fig.1| TIC 241249530 b transit measurements. a, TESS photometric 
measurements (blue), shown at their native 30-min cadence, and the best- 
fit transit model (black curve) and 30 confidence region (grey). b, Residuals 
to the best-fit model for the TESS transit detection. c, Diffuser-assisted 
ARCTIC photometric measurements (red) and the best-fit transit model 
(black curve) and 30 confidence region (grey). We show both the raw 30-s 


and RV measurements is best explained by a massive exoplanet on 
a long-period, eccentric orbit. 

To characterize the host star TIC 241249530, we separately ana- 
lysed the NEID and HARPS-N spectra using synthetic spectral fitting 
techniques and we then fit the spectral energy distribution (SED; 
see Methods). TIC 241249530 is a main-sequence star that is slightly 
hotter, larger and more massive than the Sun; the derived parameters, 
listed in Extended Data Table 1, suggest that the star is 3.2 + 0.5 Gyr old. 
The star also has a low-mass binary stellar companion, TIC 241249532, 
at a projected separation of 4.930 + 0.104”, or 1,664 +11 AU. 

We jointly fit the NEID, HPF, HARPS-N, TESS and ARCTIC measure- 
ments, accounting for perturbations to the in-transit RV signal owing 
to the Rossiter-McLaughlin effect. The transit and RV fits are shown 
in Figs. 1 and 2, respectively, and the best-fit parameters are given in 
Table 1. TIC 241249530 b is an exoplanet that is 4.98+9:1 times as mas- 
sive as Jupiter and it is on a 165.77190+9:90027-day orbit around its 
host star, with an eccentricity of 0.9412*?-00°" Our fit to the Rossiter- 
McLaughlin signal (Fig. 2b) shows that the exoplanet is orbiting inthe 
opposite direction to the projected stellar spin (A = 163.5°34) and is 
retrograde to 99.5% confidence. Few exoplanets have orbits as extreme 
as this; this orbit is more eccentric than that of any other transiting 
exoplanet, and only a handful of known planets have similarly large 
projected spin-orbit misalignments”. 

The planet that most closely resembles TIC 241249530 b is HD 
80606 b (ref. 4), which has a mass 4.1 times that of Jupiter and is also 
ona misaligned orbit with a period of 111 days and an eccentricity of 
0.93. HD 80606 bis an archetypal example of an exoplanet destined to 
become ahotJupiter with an eventual orbital period of less than 10 days. 
The eccentric orbit carries the planet close enough to its host star at 
periastron that tides raised onthe planet and star will sap energy from 
the orbit, causing it to gradually shrink and circularize. Also, simula- 
tions* show that the present orbit of HD 80606 bis consistent witha 
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cadence and binned 15-min cadence measurements, along with the 
residuals to the fitted transit signal. d, Residuals to the best-fit model for the 
ARCTIC transit detection. All brightnesses are given in parts per thousand 
(ppt). Error bars on individual data points indicate the 1o measurement 
uncertainties. 


history of von Zeipel-Lidov-Kozai (vZLK) eccentricity oscillations”? 


driven by angular momentum exchange with HD 80607, the stellar com- 
panion tothe host star. Our own simulations of the dynamical history 
and trajectory of TIC 241249530 b (see Methods) show that the orbit is 
consistent with this same type of perturber-coupled, high-eccentricity 
tidal migration. Eccentricity oscillations would have continued until the 
most recent few hundred million years, at which point general relativis- 
tic precession overtook the torque exerted by the companion, locking 


Table 1 | TIC 241249530b system parameters 


Parameter Value Description 

To 2458860.800770-.00% Time of mid-transit (BJD) 

P 165.77190+2:90927 Orbital period (days) 

e 0.9412+9:0009 Orbital eccentricity 

@ 42.32°0-40 Argument of periastron (°) 

i 85.17:9:37 Orbital inclination (°) 

K 463.374, RV semi-amplitude (ms”) 

M, 4.98+9:18 Exoplanet mass (M,) 

R, 1.19:2:94 Exoplanet radius (R,) 

A 163.5'24 Projected spin-orbit obliquity (°) 
M, 1:2712:061 Stellar mass (M,) 

R, 1.397382 Stellar radius (Ro) 

vsini, 4.607038 Projected rotational velocity (kms") 


We report the median values of the posterior distributions from our joint fit to the observed 
transits and RVs. The uncertainties represent the 68% confidence intervals (+10) for each 


parameter. 
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Fig. 2|Phase-folded RV measurements for TIC 241249530. a, RV 

measurements from NEID (blue), HPF (black), HARPS-N (yellow) and best-fit 


orbit model (black curve). b, Residuals to the RV orbit fit. c, Best-fit Rossiter- 
McLaughlin model (solid curve) and lo confidence region (grey), with the signal 


the exoplanet onan eccentric orbit that is now gradually circularizing. 
The architectures of the HD 80606 and TIC 241249530 systems lend 
support to this process asa plausible hot-Jupiter-formation mechanism. 
However, although other giant exoplanets on tidal migration tracks 
have been discovered”, including two that probably have vZLK-driven 
dynamical histories”, no previous examples have eccentricities >0.9 
and none have formation scenarios as clear as that of HD 80606 b. 
The observed occurrence rate of super-eccentric progenitors to hot 
Jupiters”** falls well short of predictions from simulations”, suggest- 
ing that giant-planet migration is dominated by other channels. With 
the discovery of TIC 241249530 b, asecond super-eccentric exoplanet 
in a hierarchical triple system has been added to the sparse sample, 
providing a new lens through which we can explore the formation of 
these planets. 

Not only do the TIC 241249530 b and HD 80606 b systems share 
similar orbital architectures but these exoplanets also have similar 
masses. The masses and eccentricities of all transiting warm Jupi- 
ters, which we define as giant planets with intermediate periods 
(10 days < P< 365 days), are shown in Fig. 3b. These two planets, 
which are the only members of the sample with super-eccentric orbits 
(e> 0.9), are also among the most massive. A correlation between 
exoplanet mass and eccentricity has been identified in previous 
works” °°, each of which found that higher-mass planets are more 
likely to have larger orbital eccentricities. We find that our narrower 
sample of transiting giant planets conforms to this known trend 
(see Methods); the eccentricity distributions of high-mass (M, > 2 M,) 
and low-mass (0.3 M, < M, < 2 M,) members of this population are 
statistically distinct (Fig. 3a). Although lower-mass planets are more 
likely to be found on low-eccentricity orbits, high-mass planets exhi- 
bit a broad, nearly flat distribution from circular to highly eccentric 
orbits. 

Although the observed mass-eccentricity correlation may be shaped 
by several processes, suchas collisional eccentricity growth” or reso- 
nant interactions with the protoplanetary disk?™”, the high masses 
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for aligned orbit shown for comparison (dashed curve). d, Residuals to the 
Rossiter-McLaughlin model fit. The dashed red box ina highlights the in- 
transit RV measurements, which are shown in b. Error bars on individual data 
points indicate the lomeasurement uncertainties. 


of TIC 241249530 b and HD 80606 b may offer a clue as to the dearth 
of super-eccentric giant planets. During the high-eccentricity phase 
of vZLK oscillations, orbital eccentricities can be driven so close to 
unity that the exoplanet will approach, or even breach, the tidal radius 
of the host star. Because the tidal radius is inversely proportional to 
the planet-star mass ratio, lower-mass planets more easily cross 
this threshold and experience tidal disruption. A relative dearth of 
low-mass, eccentric progenitors to hot Jupiters is a consistent out- 
come of simulations of high-eccentricity tidal migration following 
vZLK oscillations under an equilibrium tide assumption® ®. For planets 
susceptible to chaotic, or diffusive, dynamical tidal evolution, whereby 
oscillations excited in the planet accelerate orbital decay, this mass 
dependence is largely erased, as low-mass planets can become decou- 
pled from the perturber before being disrupted**. However, chaotic 
tides facilitate circularization on a much shorter timescale (<100 Myr); 
these planets will spend very little time with intermediate-period 
orbits****. It is possible that only the most massive eccentric giant 
planets last long enough in this period regime to be represented in 
the observed sample. 

TIC 241249530 b passes through periastron just six hours before 
each transit, presenting a unique opportunity to observe how an exo- 
planet atmosphere responds toa rapid, extreme heating event. Tem- 
poral variations in exoplanet atmospheres are best explored through 
studies of planets on eccentric orbits, for which we may see signa- 
tures of time-varying irradiation and changing pressure-temperature 
profiles, such as turbulent surface flows” and disequilibrium chem- 
istry”, depending on the heat-redistribution timescales. The atmos- 
pheres of several eccentric giant planets have been studied” ™, but the 
periastron phase has not been captured in transit for these systems. 
The orbital geometry of TIC 241249530 b will make such measure- 
ments possible for the first time. The planetary atmosphere can also 
be studied by means of emission measurements during other orbital 
phases, but the orientation precludes a secondary eclipse at a 6.30 
confidence. 
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Fig. 3 | Mass-eccentricity distribution for transiting warm Jupiters. 

a, Best-fit beta distributions for transiting warm Jupiters more massive (red) 
and less massive (blue) than 1.935 M,. The shaded regions represent the 

lo (darkest), 2oand 30 (lightest) posteriors on each fit. b, Masses and 
eccentricities of transiting warm Jupiters. The super-eccentric hot-Jupiter 
progenitors HD 80606 band TIC 241249530 bare labelled with upright and 
inverted stars, respectively. The horizontal dashed line indicates the median 
mass of the population, which is the threshold used for the fits shown ina. Our 
results are insensitive to the exact value of the threshold (see Methods); the 
median is chosen simply for visualization purposes. Error bars on individual 
data points indicate the louncertainties from the literature. 


Online content 


Any methods, additional references, Nature Portfolio reporting summa- 
ries, source data, extended data, supplementary information, acknowl- 
edgements, peer review information; details of author contributions 
and competing interests; and statements of data and code availability 
are available at https://doi.org/10.1038/s41586-024-07688-3. 


1. Dawson, R. |. & Johnson, J. A. Origins of hot Jupiters. Annu. Rev. Astron. Astrophys. 56, 
175-221 (2018). 

2. Holman, M., Touma, J. & Tremaine, S. Chaotic variations in the eccentricity of the planet 
orbiting 16 Cygni B. Nature 386, 254-256 (1997). 

3. Fabrycky, D. & Tremaine, S. Shrinking binary and planetary orbits by Kozai cycles with 
tidal friction. Astrophys. J. 669, 1298-1315 (2007). 

4. Naef, D. et al. HD 80606 b, a planet on an extremely elongated orbit. Astron. Astrophys. 
375, L27-L30 (2001). 

5. Wu, Y. & Murray, N. Planet migration and binary companions: the case of HD 80606b. 
Astrophys. J. 589, 605-614 (2003). 

6. Petrovich, C. Steady-state planet migration by the Kozai-Lidov mechanism in stellar 
binaries. Astrophys. J. 799, 27 (2015). 

7. Anderson, K. R., Storch, N. I. & Lai, D. Formation and stellar spin-orbit misalignment of hot 
Jupiters from Lidov-Kozai oscillations in stellar binaries. Mon. Not. R. Astron. Soc. 456, 
3671-3701 (2016). 

8. Muñoz, D. J., Lai, D. & Liu, B. The formation efficiency of close-in planets via Lidov-Kozai 
migration: analytic calculations. Mon. Not. R. Astron. Soc. 460, 1086-1093 (2016). 

9. Ricker, G. R. et al. Transiting Exoplanet Survey Satellite (TESS). J. Astron. Telesc. Instrum. 
Syst. 1, 014003 (2015). 


4 | Nature | www.nature.com 


10. Scott, N. J. et al. The NN-explore Exoplanet Stellar Speckle Imager: instrument description 
and preliminary results. Publ. Astron. Soc. Pac. 130, 054502 (2018). 

11. Schwab, C. et al. Design of NEID, an extreme precision Doppler spectrograph for WIYN. 
Proc. SPIE 9908, 99087H (2016). 

12. Mahadevan, S. et al. The Habitable-zone Planet Finder: a stabilized fiber-fed NIR 
spectrograph for the Hobby-Eberly Telescope. Proc. SPIE 8446, 844615 (2012). 

13. Cosentino, R. et al. Harps-N: the new planet hunter at TNG. Proc. SPIE 8446, 84461V 
(2012). 

14. Peluso, D. O. et al. The Unistellar Exoplanet Campaign: citizen science results and 
inherent education opportunities. Publ. Astron. Soc. Pac. 135, 015001 (2023). 

15. Stefánsson, G. et al. Toward space-like photometric precision from the ground with 
beam-shaping diffusers. Astrophys. J. 848, 9 (2017). 

16. Huehnerhoff, J. et al. Astrophysical Research Consortium Telescope Imaging Camera 
(ARCTIC) facility optical imager for the Apache Point Observatory 3.5m telescope. Proc. 
SPIE 9908, 99085H (2016). 

17. Albrecht, S. H., Dawson, R. I. & Winn, J. N. Stellar obliquities in exoplanetary systems. 
Publ. Astron. Soc. Pac. 134, 082001 (2022). 

18. von Zeipel, H. Sur l'application des séries de M. Lindstedt a l'étude du mouvement des 
comètes périodiques. Astron. Nachr. 183, 345 (1910). 

19. Lidov, M. L. The evolution of orbits of artificial satellites of planets under the action of 
gravitational perturbations of external bodies. Planet. Space Sci. 9, 719-759 (1962). 

20. Kozai, Y. Secular perturbations of asteroids with high inclination and eccentricity. Astron. 
J 67, 591-598 (1962). 

21. Dong, J. et al. TOI-3362b: a proto hot Jupiter undergoing high-eccentricity tidal 
migration. Astrophys. J. Lett. 920, L16 (2021). 

22. Barbieri, M. et al. HD 17156b: a transiting planet with a 21.2-day period and an eccentric 
orbit. Astron. Astrophys. 476, L13-L16 (2007). 

23. Santerne, A. et al. SOPHIE velocimetry of Kepler transit candidates. XII. KOI-1257 b: 

a highly eccentric three-month period transiting exoplanet. Astron. Astrophys. 571, A37 
(2014). 

24. Dawson, R. 1., Murray-Clay, R. A. & Johnson, J. A. The photoeccentric effect and proto-hot 
Jupiters. Ill. A paucity of proto-hot Jupiters on super-eccentric orbits. Astrophys. J. 798, 
66 (2015). 

25. Jackson, J. M. et al. Statistical analysis of the dearth of super-eccentric Jupiters in the 
Kepler sample. Astron. J 165, 82 (2023). 

26. Socrates, A. et al. Super-eccentric migrating Jupiters. Astrophys. J. 750, 106 (2012). 

27. Butler, R. P. et al. Catalog of nearby exoplanets. Astrophys. J. 646, 505-522 (2006). 

28. Ford, E.B. & Rasio, F. A. Origins of eccentric extrasolar planets: testing the planet-planet 
scattering model. Astrophys. J. 686, 621-636 (2008). 

29. Wright, J. T. et al. Ten new and updated multiplanet systems and a survey of exoplanetary 
systems. Astrophys. J. 693, 1084-1099 (2009). 

30. Frelikh, R. et al. Signatures of a planet-planet impacts phase in exoplanetary systems 
hosting giant planets. Astrophys. J. Lett. 884, L47 (2019). 

31. Papaloizou, J. C. B., Nelson, R. P. & Masset, F. Orbital eccentricity growth through 
disc-companion tidal interaction. Astron. Astrophys. 366, 263-275 (2001). 

32. Goldreich, P. & Sari, R. Eccentricity evolution for planets in gaseous disks. Astrophys. J. 
585, 1024-1037 (2003). 

33. Romanova, M. M. et al. Eccentricity growth of massive planets inside cavities of 
protoplanetary discs. Mon. Not. R. Astron. Soc. 523, 2832-2849 (2023). 

34. Vick, M., Lai, D. & Anderson, K. R. Chaotic tides in migrating gas giants: forming hot 
and transient warm Jupiters via Lidov-Kozai migration. Mon. Not. R. Astron. Soc. 484, 
5645-5668 (2019). 

35. Wu, Y. Diffusive tidal evolution for migrating hot Jupiters. Astron. J 155, 118 (2018). 

36. Rozner, M. et al. Inflated eccentric migration of evolving gas giants | - accelerated 
formation and destruction of hot and warm Jupiters. Astrophys. J. 931, 10 (2022). 

37. Langton, J. & Laughlin, G. Hydrodynamic simulations of unevenly irradiated Jovian 
planets. Astrophys. J. 674, 1106-1116 (2008). 

38. Mayorga, L. C. et al. Variable irradiation on 1D cloudless eccentric exoplanet 
atmospheres. Astrophys. J. 915, 41 (2021). 

39. Laughlin, G. et al. Rapid heating of the atmosphere of an extrasolar planet. Nature 457, 
562-564 (2009). 

40. Lewis, N. K. et al. Orbital phase variations of the eccentric giant planet HAT-P-2b. 
Astrophys. J. 766, 95 (2013). 

41. de Wit, J. et al. Direct measure of radiative and dynamical properties of an exoplanet 
atmosphere. Astrophys. J. Lett. 820, L33 (2016). 


Publisher's note Springer Nature remains neutral with regard to jurisdictional claims in 
published maps and institutional affiliations. 


Open Access This article is licensed under a Creative Commons Attribution 

a 4.0 International License, which permits use, sharing, adaptation, distribution 

and reproduction in any medium or format, as long as you give appropriate 

credit to the original author(s) and the source, provide a link to the Creative Commons licence, 
and indicate if changes were made. The images or other third party material in this article are 
included in the article's Creative Commons licence, unless indicated otherwise in a credit line 
to the material. If material is not included in the article’s Creative Commons licence and your 
intended use is not permitted by statutory regulation or exceeds the permitted use, you will 
need to obtain permission directly from the copyright holder. To view a copy of this licence, 
visit http://creativecommons.org/licenses/by/4.0/. 


© The Author(s) 2024 


'U.S. National Science Foundation National Optical-Infrared Astronomy Research Laboratory 
(NSF NOIRLab), Tucson, AZ, USA. 7Department of Astronomy and Astrophysics, The 
Pennsylvania State University, University Park, PA, USA. °Center for Exoplanets and Habitable 
Worlds, The Pennsylvania State University, University Park, PA, USA. “Department of Physics, 


Massachusetts Institute of Technology, Cambridge, MA, USA. °Kavli Institute for Astrophysics 
and Space Research, Massachusetts Institute of Technology, Cambridge, MA, USA. °Center for 
Computational Astrophysics, Flatiron Institute, New York, NY, USA. "Van Vleck Observatory, 
Astronomy Department, Wesleyan University, Middletown, CT, USA. ®Instituto de Astrofísica de 
Canarias (IAC), La Laguna, Tenerife, Spain. (Departamento de Astrofisica, Universidad 

de La Laguna (ULL), La Laguna, Tenerife, Spain. "Department of Physics and Astronomy, 
University of Pennsylvania, Philadelphia, PA, USA. "Earth and Planets Laboratory, Carnegie 
Institution for Science, Washington DC, USA. "Department of Astronomy, Indiana University 
Bloomington, Bloomington, IN, USA. “Department of Physics and Astronomy, Vanderbilt 
University, Nashville, TN, USA. “Department of Physics and Astronomy, University of New 
Mexico, Albuquerque, NM, USA. Penn State Extraterrestrial Intelligence Center, The 
Pennsylvania State University, University Park, PA, USA. School of Mathematical and Physical 
Sciences, Macquarie University, North Ryde, New South Wales, Australia. "The Macquarie 
University Astrophysics and Space Technologies Research Centre, Macquarie University, North 
Ryde, New South Wales, Australia. “Department of Astronomy and Steward Observatory, 
University of Arizona, Tucson, AZ, USA. Carl Sagan Center, SETI Institute, Mountain View, CA, 
USA. NASA Goddard Space Flight Center, Greenbelt, MD, USA. 7'Center for Planetary Systems 


Habitability, The University of Texas at Austin, Austin, TX, USA. McDonald Observatory, The 
University of Texas at Austin, Austin, TX, USA. Department of Astronomy and Astrophysics, 
University of California, Santa Cruz, Santa Cruz, CA, USA. “Jet Propulsion Laboratory, California 
Institute of Technology, Pasadena, CA, USA. “Physics Department, Hobart and William Smith 
Colleges, Geneva, NY, USA. “Department of Astronomy, Cornell University, Ithaca, NY, USA. 
7Centre for Astrophysics, University of Southern Queensland, Toowoomba, Queensland, 
Australia. Department of Astronomy and Astrophysics, Tata Institute of Fundamental Research, 
Mumbai, India. Institute of Astronomy, Faculty of Physics, Astronomy and Informatics, Nicolaus 
Copernicus University, Torun, Poland. °°Department of Physics, Lehigh University, Bethlehem, 
PA, USA. Department of Physics & Astronomy, University of California, Irvine, Irvine, CA, USA. 
“Schmidt Sciences, New York, NY, USA. Department of Astrophysical Sciences, Princeton 
University, Princeton, NJ, USA. “Anton Pannekoek Institute for Astronomy, University of 
Amsterdam, Amsterdam, The Netherlands. Department of Physics & Astronomy, University of 
Minnesota Duluth, Duluth, MN, USA. NASA Headquarters, Washington DC, USA, . Department 
of Physics, Engineering & Astronomy, Stephen F. Austin State University, Nacogdoches, TX, 
USA. Unaffiliated: Bruno Guillet, Rachel Knight, Liouba Leroux, Margaret Loose, Michael 
Primm, Masao Shimizu, Georges Simard, Stefan Will. “e-mail: arvind.gupta@noirlab.edu 


Nature | www.nature.com | 5 


Article 


Methods 


TESS photometry 

TIC 241249530 was observed with TESS? from 24 December 2019 to 
21 January 2020 (Sector 20) at a30-min cadence and from 21 Decem- 
ber 2022 to 18 January 2023 (Sector 60) at a 2-min cadence. A single 
transit-like dip (flux depth about 8 parts per thousand) was identified 
by the TESS Single Transit Planet Candidate Working Group (TSTPC 
WG) in the Sector 20 Quick Look Pipeline*?*? light curve using a box 
least-squares search. The TSTPC WG focuses on searching full-frame 
TESS light curves for isolated transit events and validating and confirm- 
ing those that are true planets, with the aim of increasing the yield of 
TESS planets with period >30 days (for example, refs. 44-47). There 
is no flux centroid motion during the transit event for TIC 241249530 
and we identify no other sources brighter than Am, = 5 in the target 
aperture. Although there is flux contamination from two nearby stars 
with 6 > Am, > 5, TIC 241249532 and TIC 241249533, both of which were 
centred onthe same pixel as TIC 241249530 in Sector 20, these are too 
faint to have been responsible for the observed change in brightness. 
No notable brightness variations were detected in the Sector 60 light 
curve. For all subsequent analysis in this work, we rely on the pre-search 
data conditioned simple aperture photometry*® °° (PDCSAP) light 
curve from the Science Processing Operations Center” (SPOC) for 
Sector 60 and the TESS-SPOC™ light curve for Sector 20 (Extended 
Data Fig. 1). 


High-contrast imaging 

To verify that the transit signature detected by TESS was indeed associ- 
ated with TIC 241249530 and not with a nearby star or binary system 
that was blended in the TESS aperture, we used the NN-EXPLORE Exo- 
planet Stellar Speckle Imager (NESSI)"° on the WIYN 3.5-m telescope 
at Kitt Peak National Observatory to conduct high-spatial-resolution 
observations of the target on 21 April 2021. A sequence of 1,000 
40-ms exposures was taken in the 832-nm and 562-nm narrow-band 
filters simultaneously with the red and blue NESSI cameras, respec- 
tively. These diffraction-limited exposures were used to reconstruct 
high-contrast images (Extended Data Fig. 2) following the steps out- 
lined in ref. 53. The achieved So contrast limits are sufficient to rule out 
the presence of faint stellar companions and background sources with 
Amag;, < 3.3 and Amag,, < 3.7 at a separation of 0.2” and Amag,., < 3.9 
and Amag,, < 4.8 at a separation of 1”. 


Ground-based photometric observations 

We used the Unistellar Network, a collaboration of citizen scientists 
using Unistellar telescopes™ in support of astronomical research, to 
observe TIC 241249530 from locations in Japan, Europe and the United 
States in search for transit signatures in March 2023. Observations were 
taken at various times from 7 to 19 March 2023, when the companion 
orbital period and transit ephemeris were still highly uncertain. After 
removing off-target and saturated frames, we calibrated the remain- 
ing images, binned them in sets of 15-30 to amplify the signal-to-noise 
ratio (S/N) and performed differential photometry®~*. No signatures of 
statistical significance were found in the Unistellar data and, based on 
our subsequent orbit fit, we confirm that none of these observations 
were taken during the transit. 

We observed TIC 241249530 again on 30 August 2023 with the Astro- 
physical Research Consortium Telescope Imaging Camera (ARCTIC)"° 
onthe ARC 3.5-m telescope at Apache Point Observatory (APO). Obser- 
vations were conducted using a beam-shaping diffuser, which creates a 
stable top-hat point spread function of the star toimprove photometric 
precision’. We used the Semrock narrow-band filter (838-876 nm) to 
avoid atmospheric absorption bands”. We began observing when the 
target rose above an air mass of 4 (altitude ~ 10°) and continued until 12° 
morning twilight, collecting a continuous 4.3-h baseline of consecutive 
30-s exposures. As the star rose above air mass approximately 1.5, about 


2.5 h after the start of the observing sequence, a transit-like decrease 
in brightness was observed. 

To reduce the ARCTIC data, a median-combined master bias image 
was constructed and subtracted from the individual science frames, 
which were flat-fielded using dome flat exposures taken at the start 
of the night. We performed differential aperture photometry on the 
reduced data using Astrolmage]* with a 17-pixel (7.7”) aperture and 
four reference stars that were carefully selected to minimize the scatter 
of the out-of-transit flux. Flux uncertainties were calculated following 
the procedures in refs. 15,59, which account for photon noise from the 
star and background, detector read noise and air-mass-dependent 
scintillation noise. We removed exposures flagged by AstrolmageJ 
for approaching the detector saturation limit, as well as exposures 
taken during intermittent cloud cover that introduced further scatter. 

The diffused point spread function of TIC 241249530 overlapped 
with that of TIC 241249532. Before initiating our ARCTIC observing 
sequence, we collected several individual exposures without the dif- 
fuser in the optical path. We used these data to calculate the relative 
brightness contributions of the two stars. TIC 241249532 contributes 
just 0.53% of the total flux in the Semrock bandpass. 


Spectroscopic observations 

We monitored the RV signal of TIC 241249530 with the NEID spectro- 
graph” onthe WIYN3.5-m telescope, collecting measurements on 40 
separate nights between 2 September 2021 and 1 March 2024. Onall but 
three of these nights, single exposures were taken, with exposure times 
ranging from 500 to 1,800 s, depending on the observing conditions. 
On the night of 30 August 2023, four consecutive 20-min exposures 
were taken simultaneously with the partial transit as observed with ARC- 
TIC, and onthe subsequent night, we secured a pair of measurements 
separated by an hour. We also obtained a sequence of 15 consecutive 
20-min exposures on the night of 12 February 2024; this sequence cov- 
ered a full transit as well as several measurements before ingress and 
after egress. We discard two spectra that were taken on nights for which 
the wavelength calibration was identified to be unreliable, leaving us 
with 56 high-quality measurements with a median S/N per extracted 
pixel of 25 at 550 nm. The raw echelle spectra were processed with ver- 
sion 1.3 of the NEID Data Reduction Pipeline (DRP; https://neid.ipac. 
caltech.edu/docs/NEID-DRP/), which produces wavelength-calibrated 
1D spectra and then calculates RVs using the cross-correlation function 
(CCF) method”. We also independently calculated the RVs from the 
calibrated 1D spectra using a modified version of the SpEctrum Radial 
Velocity AnaLyser (SERVAL) template-matching algorithm” that has 
been optimized for NEID spectra as described by ref. 63. The SERVAL 
RVs were calculated using the central 7,000 pixels of 79 orders centred 
between 4,010 and 8,400 A (order indices 20 to 100, corresponding to 
echelle orders 153 to 73). The template-matching results outperform 
the CCF-based RVs from the NEID DRP, with median single measurement 
precisions Of Ogy scavay = 6-3 MS ‘ANd Ogy pgp = 7.9 m s7, so we chose to use 
the SERVAL RVs for the analysis performed in this work. 

Further RV measurements were taken with the Habitable-zone 
Planet Finder (HPF) spectrograph”, which is on the Hobby-Eberly 
Telescope (HET)°*® at McDonald Observatory, and the HARPS-N 
spectrograph, mounted on Telescopio Nazionale Galileo (TNG) inLa 
Palma, as TIC 241249530 approached periastron in March 2023. Six 
HPF observations were made between 6 and 31 March 2023, for which 
each observation consisted of two consecutive 945-s exposures witha 
median nightly binned S/N per extracted pixel of 137 at 1,000 nm. These 
data were processed using the HxRGproc® and barycorrpy” packages 
and the RVs were calculated using a version of SERVAL that has been 
modified for HPF°*. We achieve a median RV measurement preci- 
sion of 15.0 m s”. We also observed the target five times with HARPS-N 
between 7 and 18 March 2023, with an exposure time of 3,300 s anda 
mean (min, max) S/N of 55 (37, 75). We reduced the data with the offline 
version of the HARPS-N data-reduction software through the Yabi web 


interface” installed at the Italian Center for Astronomical Archives Data 
Center. To extract the RVs, we used a G2 mask template and obtained 
a CCF width of 9.9 km s+, with an average precision of 0.1 km s”. The 
median resulting RV measurement precision is 3.4 ms. We show the 
complete RV time series from NEID, HPF and HARPS-N in Extended 
Data Fig. 3. 


Stellar characterization of TIC 241249530 

To determine the stellar atmospheric parameters of TIC 241249530, we 
analysed the out-of-transit NEID spectra collected before September 
2023 (cumulative S/N = 100 at 550 nm) using the iSpec””” Python pack- 
age to perform synthetic spectral fitting. We used the SPECTRUM radia- 
tive transfer code”, MARCS atmospheric models”, solar abundances 
from 3D hydrodynamic models” and the sixth version of the Gaia ESO 
survey (GES) atomic line list”®. The microturbulence velocity was treated 
as a free parameter to allow for flexibility in accounting for small-scale 
motions inthe stellar atmosphere. Macroturbulence was determined 
using an empirical relation, making use of established correlations with 
other stellar properties”. To streamline the fitting, we restricted the 
analysis to specific spectral regions from 480 to 680 nm, encompassing 
the wing segments of the Ha, HB and Mg I triplet lines, which are sensi- 
tive to 7,,,and logg, and the Fe I and Fe II lines, which provide precise 
constraints on [Fe/H] and vsini,. We minimize the difference between 
the synthetic and input spectra by applying the nonlinear least-squares 
Levenberg-Marquardt fitting algorithm, using constraints from the 
aforementioned models and line lists. 

The HARPS-N spectra were independently analysed with BACCHUS®, 
using MARCS atmospheric models, the GES atomic line list and the TUR- 
BOSPECTRUM radiative transfer code”*®®. For our fit, we constrained Ter 
by requiring Fe I line abundances to be uncorrelated with their respec- 
tive excitation potentials in the synthetic spectrum and we constrained 
logg by requiring ionization balance for the Fe land Fe II lines. We also 
required the Fe I line abundances to be uncorrelated with their equiva- 
lent widths and the stellar metallicity ([Fe/H]) was calculated as the 
average of these abundances. The projected rotational velocity was 
estimated by fitting the broadening of the Fe I lines, accounting for 
the best-fit microturbulence and assuming the same macroturbulence 
contribution as in the iSpec analysis. The stellar parameters derived 
from the NEID and HARPS-N spectra are largely in good agreement 
(<lo). Discrepancies between the [Fe/H] values and vsini, values at the 
1.20 level probably result from differences between the fitted microtur- 
bulence, which is known to exhibit small variations for different fitting 
methods”. We adopt the iSpec Ty, logg, [Fe/H] and vsini, for the rest 
of the analysis in this work. 

We performed an analysis of the broadband SED of TIC 241249530 
together with the Gaia DR3 parallax following the procedures described 
in refs. 81-83. We use JHK; magnitudes from 2MASS*, W1-W3 magni- 
tudes from WISE®, GgpGgp magnitudes from Gaia®, BVgri magnitudes 
from APASS® and the NUV magnitude from GALEX®*. We also used 
the Gaia spectrophotometry spanning 0.4-1.0 um. Altogether, the 
available photometry spans the full stellar SED over the wavelength 
range 0.2-10.0 um. We fit the SED using PHOENIX stellar atmosphere 
models®’, with the effective temperature, surface gravity and metal- 
licity set to the spectroscopically determined values. The remaining 
free parameter is the extinction (A,), which we limited to the maximum 
line-of-sight value of A, = 0.44 mag from galactic dust maps”’. The 
resulting fit is shown in Extended Data Fig. 4. Integrating the unred- 
dened model SED yields the bolometric flux at Earth, F,,,, = 7.19 + 0. 
20 x 10™ erg scm”. Taking the Fpa and Tese together with the Gaia 
parallax, we calculate the stellar radius to be R, = 1.404 + 0.028 Rg. 
Also, the stellar mass is inferred using empirical relations”, giving 
M, =1.24 + 0.07 M,, and we estimate the age to be 3.2 + 0.5 Gyr by fit- 
ting the evolutionary state with the Yonsei-Yale isochrone models”. 
Our reported 0.5-Gyr uncertainty accounts for the uncertainties on 
each of the inputs to the isochrone fit: effective temperature, surface 


gravity, metallicity and stellar mass. However, this does not account 
for systematic uncertainties arising from our choice of stellar models, 
which can be on the order of 1 Gyr. 

The best-fit extinction for our SED model is Ay = 0.31 + 0.02. This 
large value is supported by a clear detection of interstellar absorption 
in the Na D doublet and the KI 770 nm lines in the NEID spectra. Both 
spectroscopic analyses yield Tę values that are substantially hotter 
than the literature value from Gaia DR3 spectrophotometric analysis”, 
which is consistent with the effect of reddening from dust along the 
line of sight to the star. 

Using the projected rotational velocity and stellar radius, we place 
an upper limit on the rotation period of 16.9738 days. We attempt to 
make amore precise measurement of the rotation period to determine 
the stellar inclination, but the existing data are insufficient. An analy- 
sis of the TESS light curves using the TESS Systematics-Insensitive 
Periodogram package” shows no notable photometric modulation on 
timescales shorter than the length of each individual sector. We also 
examine archival photometry of the star from the WASP survey”. These 
data consist of 2,178 measurements on 33 nights, with two isolated 
epochs in April 2006 and March 2007, and the remaining data covering 
October 2007 to March 2008. In spite of the substantially longer base- 
line than TESS, a _Lomb-Scargle periodogram analysis of the WASP 
measurements reveals no notable peaks besides the half-day, one-day 
and two-day sampling aliases. The lack of photometric modulation is 
reflected in the spectroscopic data as well; we do not detect periodic 
variation in the activity-sensitive Ca II H & K, Nal or Ha spectral lines 
as measured by the NEID DRP. Also, there is no emission in the Ca II H 
& K line cores in the NEID and HARPS-N spectra, suggesting that the 
star is chromospherically quiet. 


Stellar characterization of TIC 241249532 

TIC 241249530 shares a common parallax and proper motion with 
TIC 241249532 as measured by Gaia, and the two stars are separated 
on the sky by 4.930 + 0.104” (ref. 86). The probability of a chance 
alignment between TIC 241249530 and TIC 241249532 is R = 9.73 x 10° 
(ref. 95), suggesting that the pair is indeed gravitationally bound. Gaia’s 
photometric measurements of TIC 241249532 place it firmly along 
the main sequence. We do not perform an independent SED analy- 
sis on this star but instead estimate its mass using empirical mass- 
luminosity relations’®”. We calculate the mass to be 0.453 + 0.012 
times that of the Sun based on the 2MASS K,-band magnitude and 
0.400 + 0.016 times that of the Sun based on the Gaia G,p»-band mag- 
nitude. The stellar mass, coordinates and broadband photometry are 
given in Extended Data Table 1. 

Onthe basis of the weighted mean of the Gaia parallax measurements 
for the system, the on-sky separation corresponds to a projected physi- 
cal separation of 1,664.00 + 10.85 AU. The relative motions of these two 
stars are not constrained well enough by Gaia to meaningfully estimate 
an orbital solution. However, as an effort to quantify the dynamical 
impact of TIC 241249532 in our analysis, we simulated 10 million orbits 
sampled randomly in phase, uniformly in cosi and thermally (f(e) = 2e) 
in eccentricity. We determine the orbital period of the system to be 
>10,000 years, witha peak in frequency at 35,000 years. Long-period 
stellar companions suchas this can directly bias RV analyses of exoplan- 
ets inthe form ofa linear RV slope. However, our simulations show that 
TIC 241249532 probably induces a linear trend in the observed RVs of 
TIC 241249530 at the level of just 1 cm s” year", with 99% of our orbits 
returning slopes <30 cm s” year”. The amplitude of this signal is small 
compared with the km-s “level variations induced by TIC 241249530 b 
and we therefore do not include it as an extra body in our joint fitting. 


Joint transit + RV analysis 

We use the exoplanet software package” to fit a transit model anda 
Keplerian orbit with Rossiter-McLaughlin perturbations to the 
observed photometric and RV signals for TIC 241249530. The exoplanet 
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package relies on starry’ and the underlying analytic models from 


ref. 101 to fit the transits, and the orbital parameter posteriors are sam- 
pled using the PyYMC3 Hamiltonian Monte Carlo package". The orbit 
model consists of a full Keplerian with tight Gaussian priors on the 
orbital period, P, and time of conjunction, To, broad uniform priors on 
the exoplanet mass, M,, transit impact parameter, b, and transit depth, 
6, and Gaussian priors on the stellar mass, M,, and radius, R,. We repa- 
rameterize the eccentricity, e, and argument of periastron, œ, as /esinw 
and ./ecosw and we sample these on the unit disk. We do not impose 
an extra eccentricity prior, as the global warm-Jupiter eccentricity 
distribution is not well constrained’. Separate quadratic 
limb-darkening coefficients, reparameterized as qand q, asin ref. 104, 
are used for each instrument for which in-transit observations were 
taken (TESS, ARCTIC, NEID). For the RV data, we fit individual zero-point 
offset terms (ø) and jitter terms (y) for each instrument, splitting the 
NEID RVs into two separate datasets before and after the instrument 
restart. Dilution terms are included for both TESS and ARCTIC, as both 
transit measurements suffered from flux contamination. The TESS 
data products already account for dilution, but previous works have 
demonstrated that these results are susceptible to overcorrection”™, 
so we allow the TESS dilution term to float uniformly from 0.1to 1.2. 
For ARCTIC, we fix the dilution to be 0.9947 based on the out-of-transit 
data for which the target was well resolved from its companion. To 
model the Rossiter-McLaughlin signal, we adopt the formalism of 
ref. 107 along with their prior distributions for the Gaussian line disper- 
sion parameter, p, the Lorentzian line dispersion parameter, y, and 
macroturbulence, ĝ. We place a Gaussian prior onthe projected stellar 
rotational velocity, vsini,, and a uniform prior on the projected spin- 
orbit misalignment, A. The prior distributions and posterior results for 
all of these fit parameters, as well as for some derived values, are given 
in Extended Data Table 2. 


Stellar obliquity 
The stellar obliquity, y, is related to the projected obliquity, A, by 


cosy = sini,sinicosA + cosi,cosi. 


Here i, is the inclination of the stellar spin axis and iis the inclination 
of the exoplanet orbit. We cannot directly calculate y because the stel- 
lar inclination is not known. Instead, assuming that the stellar inclina- 
tion is drawn from an isotropic distribution, uniform in cosi,, we use 
the above equation to determine the possible values of y and their 
relative probabilities. For our derived posteriors onA and i, we find that 
the orbit is indeed retrograde (that is, y > 90°) at 99.5% confidence in 
this scenario, and we calculate the obliquity to be y = 14173). This value 
is consistent with expectations for vZLK-driven migration; simulations 
show that the final obliquity can be as large as 180° for systems such 
as this’. This is not definitive proof of the formation history, however, 
as retrograde orbits suchas this can also be produced through planet- 
planet interactions”. Regardless, we warn that our result is strongly 
dependent on the naive assumption of an isotropic stellar inclination 
distribution, which is not always valid!°™, 


Dynamical history—analytic constraints 
The high eccentricity and tight orbit of TIC 241249530 b and the pres- 
ence of the distant stellar binary companion indicate a likely history 
of high-eccentricity migration driven by vZLK oscillations and tidal 
dissipation. To determine how this formation channel could have deliv- 
ered the exoplanet to its current orbit, we first identify a set of initial 
conditions consistent with the present-day architecture of the system. 
We workin the context of the secular approximation for the evolution 
of hierarchical triple configurations. 

The planet is at present close enough to the primary star such 
that short-range forces—general relativity (GR), tides and rotational 
distortions—have quenched any vZLK oscillations driven by the 


companion. We calculate the semimajor axis at which this quenching 
occurs by assuming that GR dominates the short-range forces and 
examining the ratio between the timescale of GR precession of the 
inner orbit and the timescale for vZLK oscillations. To leading order 
(quadrupole limit), this ratio is” 


tor = at (1- e7)M,c? 
tuad 303 (1- e”)*7G(M, +M)? 


in which M, is the mass of the primary star; M,, a and e are the mass, 
semimajor axis and eccentricity of the planet; and M,, a, and e, are the 
mass, semimajor axis and eccentricity of the binary companion. For the 
exoplanet and primary-star parameters, we adopt the median posteri- 
ors from our joint fit. For the binary companion, we set M, = 0.453 M,, 
a, = 1,664 AU and e, = 0.5. If the planet started with a low initial eccen- 
tricity ofe = 0.1, vZLK oscillations would have started only if the initial 
semimajor axis of the planet was a, > 4.2 AU, for which this value is cal- 
culated by setting t¢p/tquaa = 1. 

We can now constrain the initial eccentricity by requiring that the 
periastron distance of the first vZLK oscillation was sufficiently small 
to trigger efficient tidal dissipation. In particular, in the quadrupole 
limit, the quantity aç = a (1 - e2,,,) is approximately conserved through- 
out the tidal migration, as the orbital angular momentum is conserved 
both during episodes of maximum eccentricity, as well as after vVZLK 
oscillations have been quenched. Here e,,,, indicates the maximum 
eccentricity reached during a vZLK oscillation and a; is equal to the 
final semimajor axis once the orbit has fully circularized. If a, is taken 
to be conserved, the maximum eccentricity of the initial vZLK oscilla- 
tion must have been e€; max > 0.9947. 

Exciting an eccentricity this high on the initial vZLK oscillation must 
have required a substantial initial inclination, /,, between the orbit of 
the planet and that of the binary companion. We derive a lower bound 
on/, using the following equation from ref. 112: 


1 9 Chicas j 5 2 
Ecr| > - 1) => 3 Vo, — 3 COSI |. 
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Here j nin = y1 e? max and we have assumed GR perturbations to be 
dominant over those from tidal and rotational distortion. The dimen- 
sionless quantity €,, measures the ‘strength’ of perturbations from GR 


relative to those of the binary companion and it is defined as 


Ecr=3G(M, + M,)?a3 (L- e3)°7/(Mya‘c?). 


Extended Data Fig. 5 shows the required initial inclination between 
the planetary and binary orbital planes and the resulting initial maxi- 
mum eccentricity with respect to the initial semimajor axis of the 
planet. Although vZLK oscillations are present when a; > 4.2 AU, not 
all values above this threshold yield defined values for /, because the 
short-range forces are too strong for the planet to reach the required 
high initial eccentricity unless the initial semimajor axis exceeds 
a,>7.0 AU. The maximum eccentricity of the initial VZLK cycle must 
have been e€; max > 0.9947 to generate the present-day semimajor axis 
and eccentricity. Attaining a maximum eccentricity this large is only 
possible with a nearly polar initial inclination between the orbit planes 
of the planet and the binary companion. We find that the initial incli- 
nation /, > 86.8° for a, > 7.0 AU. Altogether, these results indicate that 
it is possible to reach the present-day parameters of the system if the 
planet started beyond a, > 7.0 AU and the binary companion started on 
an orbit nearly perpendicular to that of the planet. 


Dynamical history—simulations 
We now use our derived constraints onthe initial orbital conditions to 
explore the planetary orbital evolution through numerical simulations. 


We conduct integrations of the secular equations of motion for 
TIC 241429530 through KozaiPy, a publicly available software package 
that simulates hierarchical three-body systems (https://github.com/ 
djmunoz/kozaipy). The equations of motion are provided in ref. 3. We 
adopt initial values of a; = 10 au and e; = 0.1, considering the minimum 
semimajor axis necessary for vZLK oscillations to be present. We con- 
sider perturbations to the octupole order and also account for tidal 
evolution inthe constant-time-lag model of equilibrium tides™. Tidal 
parameters are adjusted so that the system reaches its present-day 
orbital parameters at an age of 3 Gyr, approximately equal to the derived 
age of the system. Specifically, the Love number of the planet is set to 
k, = 0.25 and its viscous timescale is set to t, = 0.01 days. 

The simulation results are presented in Extended Data Fig. 6. They 
indicate the presence of vZLK oscillations that trigger periods of very 
large eccentricities. At the times that the periapse distance is mini- 
mized, tidal dissipation is strong and the semimajor axis shrinks. Even- 
tually, the semimajor axis becomes small enough that the vZLK 
oscillations are suppressed owing to short-range forces and the planet 
decouples fromthe binary companion. After the vZLK oscillations are 
quenched, the mutual inclination is approximately conserved and the 
eccentricity of the planet slowly damps owing to continued tidal dis- 
sipation. We observe that there is an instant in time at which the eccen- 
tricity and semimajor axis of the planet are very close to the present-day 
values. We also note that the value of a(1 - e2,,,) is conserved during 
episodes of maximum eccentricity of each individual vZLK cycle, rang- 
ing within only a few percent of the average value of a(1 - e2,,,). Accord- 
ing to this simulation, continued dissipation will cause the planet to 
reachacircular orbit in about a billion years. Altogether, this simulation 
provides a plausible proof of concept of the system’s dynamical history 
of coupled vZLK oscillations and tidal migration. We suggest future 
work onthe system to explore the role that dynamical tides might have 
played in its formation. 


Modelling the transiting-warm-Jupiter eccentricity distribution 
To explore the relation between exoplanet mass and eccentricity for 
warm Jupiters, or intermediate-period giant planets, we start with the 
sample ofall transiting exoplanets with masses between 0.3 and 15 times 
that of Jupiter and orbital periods between 10 and 365 days. For each sys- 
tem inthis sample, we adopt the most up-to-date mass and eccentricity 
constraints for which the eccentricity was fit as a free parameter when 
solving for the orbit. We discard four planets for which all literature 
solutions assumed a circular orbit with the eccentricity fixed at O. All 
of these planets are less than 1.3 Jupiter masses. We also remove two 
planets in P-type circumbinary orbits, as the dynamical environments 
of these systems are expected to differ from those of planets orbiting 
single stars™"ć, Our sample differs from that analysed in ref. 41, which 
draws fromthe RV planet sample and thus uses projected planet mass 
(M, sini) instead of true mass. By restricting our analysis to transiting 
exoplanets, we ensure that the measured masses are not degenerate 
with orbital inclination. This approach also mitigates the susceptibil- 
ity of our results to detection biases, as the completeness fractions of 
transit surveys should be largely insensitive to exoplanet mass in the 
Jupiter-sized-planet regime. 

The median mass of our sample is 1.935 Jupiter masses. We divide 
the sample into two groups of equal size, placing planets more massive 
than the median into one group and planets less massive than the 
median into the other group. The population-level eccentricity distri- 
bution of each group is then modelled in PyMC* using a hierarchical 
Bayesian framework. For our model, we adopt a beta distribution with 
two hyperparameters, 0 = {u, x}, in which y describes the mean of the 
distribution and 1/x describes its variance. These hyperparameters 
represent a reparameterization of the standard beta distribution 
parameters a and $, in which a = uk and B= (1 - w)x. Abeta distribution 
is chosen for its flexibility in shape and because it is naturally bounded 
between 0 and 1. We adopt a uniform hyperprior for u ~ U(O, 1) anda 


log-normal hyperprior for log x ~ N(O, 3). These choices reduce the 
impact of hyperprior choices onthe inference results, especially when 
the sample size is small, as is the case here”. The best-fit distributions 


are shownin Fig. 3 and the resulting hyperparametersare ji, = 0.187005, 


l0gKiow = 1.237935» Hien = 0-44 76.08 and l0gKhign = 0.917934. The mean 


values, y, of the eccentricity distributions of low-mass and high-mass 
transiting warm Jupiters differ by 4.20. To assess the robustness of this 
result, we repeated the process for mass cutoffs between 1Jupiter mass 
and 2.7 Jupiter masses. These bounds were chosen such that the size 
ratio of the two groups does not exceed 2:1. At each cutoff, we ran 1,000 
trials, drawing the planet masses from asymmetric Gaussian distribu- 
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uncertainties. For all mass-cutoff values over this range, the mean val- 
ues of the two eccentricity distributions differ by 3-So. 
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Extended Data Fig. 1| TESS light curves for TIC 241249530. a, Sector 20 TESS-SPOC data (30-min cadence). The vertical dashed line marks the best-fit transit 
midpoint. b, Sector 60 SPOC data (2-min cadence). Error bars on individual points represent 1o measurement uncertainties. 
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Extended Data Fig. 2 | Reconstructed NESSI speckle images and 50 contrast 
curves for TIC 241249530. Observations were taken simultaneously at 562 nm 
with the blue camera (upper-left inset image) and at 832 nm with the red camera 
(upper-right inset image). The contrast curves indicate the limiting magnitude 
difference at which bound or background sources could be detected for 
separations between 0.2” and1.2”. 
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Extended Data Fig. 3 | RV time series for TIC 241249530. a, RV measurements confidence intervals for the fit. The horizontal axis is in units of days relative to 
from NEID (blue), HPF (black) and HARPS-N (gold) and best-fit orbit model BarycentricJulian Date (BJD) 2457000. Error bars on individual points 

(black curve). b, Residuals to the RV orbit fit. Vertical red lines in both panels represent lomeasurement uncertainties. 

mark the predicted transit times. The grey-shaded region bounds the 3a 


—-9.0 


y 
-9.5 z 
"m 
| 
wn 
a -10.0 5 
E£ f = 
9 uw 
7 Z 
wu D 
D -10.5 Ò 
ace 
e 
% -11.0 
D 
fe) 
5-115 
x 
-12.0 
-12.5 
1071 10° 107 
Wavelength (A) [um] 
Extended Data Fig. 4| SED for TIC 241249530. Red symbols represent the the effective width of the passband. Blue symbols are the model fluxes from the 
observed photometric measurements, for which the vertical error bars best-fit PHOENIX atmosphere model (black). The Gaia spectrophotometry is 


represent the lomeasurement uncertainties and the horizontal bars represent represented asa grey swathe (see also inset plot). 
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the present-day orbit owing to strong short-range forces (primarily GR, with 

further contributions from tides and rotational distortions). The blue region 
(a; > 7 AU) indicates the presence of vZLK oscillations that can drive the planet 


to reach present-day parameters. 


Extended Data Fig. 5 | Constraints on the parameters of the first vVZLK 
oscillation for TIC 241249530 b. a, Constraints on initial inclination. 

b, Constrains on initial maximum eccentricity. The red region (0 AU < a, < 4.2 AU) 
indicates the absence of vZLK oscillations. The orange region (4.2 AU < a; < 7 AU) 
indicates the presence of vZLK oscillations but with insufficient e; max tO reach 
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Extended Data Fig. 6 | Simulated evolution of the orbit of TIC 241249530b 
resulting from high-eccentricity migration driven by vZLK oscillations. 

a, Evolution of the eccentricity of the planetary orbit over time. b, Evolution 

of the mutual inclination between the planet and binary orbits over time. 

c, Evolution of the semimajor axis, a, and periastron separation, a(1 - e), of the 
planetary orbit over time. The vertical red lines at 3.2 Gyr mark the age at which 
the orbit reaches the present-day conditions (e = 0.94, a= 0.64 AU). We adopt 
a,=10 AU, e, = 0.1 for the initial orbital parameters for illustrative purposes. 
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Extended Data Table 1| Stellar parameters for TIC 241249530 and TIC 241249532 


Parameter 


Identifiers: 
TIC 
Gaia DR3 
2MASS 


Coordinates and Parallax: 
R.A. 
Dec 
Parallax 
Proper Motion R.A. 
Proper Motion Dec 


Broadband photometry: 


Derived Stellar Parameters: 


Teff 


log g 
[Fe/H] 


vsin i, 


Age 
Av 


Data from refs. 84-87,119. 


Value 


241249530 


995754258812449024 


J062114114+5317423 


06:21:14.02 
+53:17:42.32 
2.96 + 0.02 
7.620 + 0.018 
-15.717 + 0.015 


12.226 + 0.413 
11.689 + 0.039 
11.1328 + 0.0079 
11.5593 + 0.0003 
11.8835 + 0.0007 
11.0698 + 0.0005 
10.538 

10.266 + 0.024 
10.192 + 0.021 
10.115 + 0.023 
10.147 + 0.021 
10.176 + 0.065 
8.993 


6166 + 43 
6150 + 50 
4.23 + 0.06 
4.25 + 0.10 
0.09 + 0.03 
0.20 + 0.09 
4.2 4 0.8 

5.5 +0.7 
1.24 + 0.07 
1.404 + 0.028 
2.56 + 0.08 
7.19 + 0.20 x 10710 
3.2 + 0.5 
0.31 + 0.02 


241249532 


995754254517548032 


J06211389+5317380 


06:21:13.76 
+53:17:37.85 
3.17 + 0.09 
7.766 + 0.107 
-15.799 + 0.092 


16.5809 + 0.0314 
17.4606 + 0.0029 
18.046 + 0.090 
16.054 + 0.026 
12.017 

13.741 + 0.064 
13.506 + 0.071 


0.453 + 0.012 


Unit 


log (cm s~?) 
log (cm s~?) 
dex 
dex 
km s71 


erg s™1 cm2 
Gyr 
mag 


85 
85 
85 
85 


iSpec 
BACCHUS 
iSpec 
BACCHUS 
iSpec 
BACCHUS 
iSpec 
BACCHUS 


Extended Data Table 2 | Joint-fit priors and posteriors 


Parameter Prior Posterior Unit 
Stellar Parameters: 
M. N (1.24, 0.07) 12710061 Mo 
R, N (1.404, 0.028) 1.397+9025 Ro 
vsin i, N (4.2, 0.8) 4.607528 km s“1 
Ç U(2.0, 6.5) 3.671, = kms 
B u(0.8, 1.2 1.017573 kms 
Y U(2.5, 4.5) 3.3773 60 km s~! 
q1, TESS u(0.0, 1.0) 0.357973 
q2TESS U(0.0, 1.0) 0.177932 
q1,ARCTIC u(0.0, 1.0) 0.477538 
2,ARCTIC U(0.0, 1.0) 0.141032 
q1,NEID u(0.0, 1.0) 0.487533 


2,NEID 026 


U(0.0, 1.0) 


Orbital and Planetary Parameters: 


—0.01+?32 


To N(2458860.8, 0.01)  2458860.800770-0018 BJD 
P N(165.77, 0.1) 165.771902 000025 days 
Ve cos w U(-1.0, 1.0) 0.717470 pone 

Ve sin w U(-1.0, 1.0) 0.65317 0.0050 

e derived 0.941279 9908 

w derived 42.321040 deg 

b U(0.0, 1.0) 0.58170 -528 

a derived 0.64179 079 AU 

A U(-180, 180) 163.5734 deg 

RY R, u(0.05, 0.15) 0.087373 3059 

Mp U(0.5, 30) 4.987076 Ms 

Rp derived 118649037 Ry 

Instrument Parameters: 

TESS Dilution u(0.1, 1.2) 1.01079 06? 

Ophot, TESS £(0.74, 10) 0.004+5:539 ppt 

Ophot, ARCTIC £(4.18, 10) 1.51745 ppt 

ORV.NEID-pre U(0, 25) 3.4128 ms 
ORV,NEID-post U(0, 25) 4.6115 m s™ 
ORV,HPF u(0, 25) 14.6134 ms 
ORV,HARPS-N U(0, 25) 3.313? m s™ 
YRVNEID-pre N (0, 1000) 57.2476 ms“ 
YRV,NEID-post N (0, 1000) —70.1733 ms 
YRVHPF (0, 1000) —18.973° m s™ 
YRVHARPS-N N (0, 1000) —170.7t5} ms” 


We report the median values of the posterior distributions from our joint fit. Uncertainties 
represent the 68% confidence intervals (+10) for each parameter. Limb darkening is sampled 
using the parameterization given by ref.104. For À, we use the custom angle distribution from 
the PyMC3 extras extension of the exoplanet package® to avoid discontinuities at +n. 


